classdef KrusellRobotsEach < handle
% KrusellRobotsEach implements the production function and its derivatives
% See Mathematica Notebook :code:`CESRobotsSkillPrices.nb` Section
% **Combination of Krusell and all labor substitutable by robots
% (KrusellRobotsEach) with labor scaled**. delta in the Mathematica
% file corresponds to scale_L here
%
% in the paper, different symbols are used. For the correspondence
% see tab_cal_prod_fun.R
    
    properties
        alpha
        rho
        sigma
        kappa_M
        kappa_R
        kappa_C
        nu
        lambda
        mu
        eta_M
        eta_R
        eta_C
        upsilon
        A
        gamma
        parse_inputs
        KrusellRobotsEach_use_properties
        scale_L
    end
    
    methods
        function self = KrusellRobotsEach(varargin)
            p = inputParser;
            p.KeepUnmatched = true;
            p.addParameter('alpha', [], @(x)  x>0)
            p.addParameter('rho', [], @(x) x>0)
            p.addParameter('sigma', [], @(x) x>0)
            p.addParameter('kappa_M', [], @(x) x>0)
            p.addParameter('kappa_R', [], @(x) x>0)
            p.addParameter('kappa_C', [], @(x) x>0)
            p.addParameter('nu', [], @(x) x>0)
            p.addParameter('lambda', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('mu', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('eta_M', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('eta_R', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('eta_C', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('upsilon', [], @(x) ((x>=0) && (x<=1)))
            p.addParameter('A', [], @(x) x>0)
            p.addParameter('gamma', [], @(x) x>0)
            p.addParameter('parse_inputs', true, @islogical)
            p.addParameter('KrusellRobotsEach_use_properties', false, @islogical)
            p.addParameter('scale_L', [], @isnumeric)
            
            p.parse(varargin{:});
        
            self.alpha = p.Results.alpha;
            self.rho = p.Results.rho;
            self.sigma = p.Results.sigma;
            self.kappa_M = p.Results.kappa_M;
            self.kappa_R = p.Results.kappa_R;
            self.kappa_C = p.Results.kappa_C;
            self.nu = p.Results.nu;
            self.lambda = p.Results.lambda;
            self.mu = p.Results.mu;
            self.eta_M = p.Results.eta_M;
            self.eta_R = p.Results.eta_R;
            self.eta_C = p.Results.eta_C;
            self.upsilon = p.Results.upsilon;
            self.A = p.Results.A;
            self.gamma = p.Results.gamma;
            self.scale_L = p.Results.scale_L;
            self.parse_inputs = p.Results.parse_inputs;
            self.KrusellRobotsEach_use_properties = ...
                p.Results.KrusellRobotsEach_use_properties;

            if self.KrusellRobotsEach_use_properties
                fprintf('setting parse_inputs to true\n')
                self.parse_inputs = true;
            end
        end
        
        
        function out = Y(self,inp,varargin)

            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_S.^par.alpha.*par.A.*((-1).*((-1)+par.mu).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^(((-1)+ ...
  par.rho).^(-1).*par.rho.*((-1)+par.sigma).*par.sigma.^(-1))+ ...
  par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((( ...
  -1)+par.sigma).*par.sigma.^(-1))).^((-1).*((-1)+par.alpha).*((-1)+ ...
  par.sigma).^(-1).*par.sigma);
        end
        
        function out = Y_M(self,inp,varargin)
            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_S.^par.alpha.*par.A.*(1+(-1).*par.alpha).*par.eta_M.* ...
  par.mu.*par.upsilon.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1) ...
  .*par.kappa_M.^(-1)).*((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1)).*(((-1).*inp.K_B_M.^(((-1)+ ...
  par.kappa_M).*par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.* ...
  par.eta_M.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1).* ...
  par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(( ...
  -1).*par.nu.^(-1)).*(par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+ ...
  par.kappa_M).*par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.* ...
  par.eta_M.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1).* ...
  par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^((( ...
  -1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).* ...
  inp.K_B_R.^(((-1)+par.kappa_R).*par.kappa_R.^(-1)).*((-1)+ ...
  par.eta_R)+inp.L_R.*par.eta_R.*par.scale_L.*(inp.L_R.* ...
  par.scale_L).^((-1).*par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^( ...
  -1).*par.kappa_R)).^(((-1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu) ...
  .^(-1)).*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((-1) ...
  .*par.sigma.^(-1)).*((-1).*((-1)+par.mu).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^(((-1)+ ...
  par.rho).^(-1).*par.rho.*((-1)+par.sigma).*par.sigma.^(-1))+ ...
  par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((( ...
  -1)+par.sigma).*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+( ...
  -1).*par.alpha.*par.sigma));

        end
        
        
        function out = Y_R(self,inp,varargin)

            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_S.^par.alpha.*par.A.*((-1)+par.alpha).*par.eta_R.*par.mu.*(( ...
  -1)+par.upsilon).*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1)).*((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1)).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(( ...
  -1).*par.nu.^(-1)).*(par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+ ...
  par.kappa_M).*par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.* ...
  par.eta_M.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1).* ...
  par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^((( ...
  -1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).* ...
  inp.K_B_R.^(((-1)+par.kappa_R).*par.kappa_R.^(-1)).*((-1)+ ...
  par.eta_R)+inp.L_R.*par.eta_R.*par.scale_L.*(inp.L_R.* ...
  par.scale_L).^((-1).*par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^( ...
  -1).*par.kappa_R)).^(((-1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu) ...
  .^(-1)).*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((-1) ...
  .*par.sigma.^(-1)).*((-1).*((-1)+par.mu).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^(((-1)+ ...
  par.rho).^(-1).*par.rho.*((-1)+par.sigma).*par.sigma.^(-1))+ ...
  par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((( ...
  -1)+par.sigma).*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+( ...
  -1).*par.alpha.*par.sigma));
        end
        
        
        function out = Y_C(self,inp,varargin)
            
            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = (-1).*inp.K_S.^par.alpha.*par.A.*((-1)+par.alpha).*par.eta_C.*(( ...
  -1)+par.lambda).*((-1)+par.mu).*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1)).*((-1).*inp.K_B_C.^(((-1) ...
  +par.kappa_C).*par.kappa_C.^(-1)).*((-1)+par.eta_C)+inp.L_C.* ...
  par.eta_C.*par.scale_L.*(inp.L_C.*par.scale_L).^((-1).* ...
  par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^(-1)).*(((-1).* ...
  inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^((-1).*par.rho.^(-1)).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^((par.rho+( ...
  -1).*par.sigma).*(par.sigma+(-1).*par.rho.*par.sigma).^(-1)).*(( ...
  -1).*((-1)+par.mu).*(inp.K_E.^(((-1)+par.rho).*par.rho.^(-1)).* ...
  par.lambda+(-1).*((-1)+par.lambda).*(((-1).*inp.K_B_C.^(((-1)+ ...
  par.kappa_C).*par.kappa_C.^(-1)).*((-1)+par.eta_C)+inp.L_C.* ...
  par.eta_C.*par.scale_L.*(inp.L_C.*par.scale_L).^((-1).* ...
  par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^(-1).*par.kappa_C)).^((( ...
  -1)+par.rho).*par.rho.^(-1))).^(((-1)+par.rho).^(-1).*par.rho.*(( ...
  -1)+par.sigma).*par.sigma.^(-1))+par.mu.*((par.upsilon.*(((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^(((-1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+ ...
  par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(((-1)+par.nu).* ...
  par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^(((-1)+par.sigma) ...
  .*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+(-1).*par.alpha.* ...
  par.sigma));
        end
        
        function out = Y_B_M(self,inp,varargin)
            
            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_B_M.^((-1).*par.kappa_M.^(-1)).*inp.K_S.^par.alpha.*par.A.*( ...
  (-1)+par.alpha).*((-1)+par.eta_M).*par.mu.*par.upsilon.*((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1)).*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).* ...
  ((-1)+par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^((-1).*par.nu.^(-1)).*(par.upsilon.*(((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^(((-1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+ ...
  par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(((-1)+par.nu).* ...
  par.nu.^(-1))).^(((-1)+par.nu).^(-1)).*((par.upsilon.*(((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^(((-1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+ ...
  par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(((-1)+par.nu).* ...
  par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((-1).* ...
  par.sigma.^(-1)).*((-1).*((-1)+par.mu).*(inp.K_E.^(((-1)+par.rho) ...
  .*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*(((-1).* ...
  inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^(((-1)+ ...
  par.rho).^(-1).*par.rho.*((-1)+par.sigma).*par.sigma.^(-1))+ ...
  par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((( ...
  -1)+par.sigma).*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+( ...
  -1).*par.alpha.*par.sigma));
        end
        
        function out = Y_B_R(self,inp,varargin)
            
            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
             out = (-1).*inp.K_B_R.^((-1).*par.kappa_R.^(-1)).*inp.K_S.^par.alpha.* ...
  par.A.*((-1)+par.alpha).*((-1)+par.eta_R).*par.mu.*((-1)+ ...
  par.upsilon).*((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1)).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(( ...
  -1).*par.nu.^(-1)).*(par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+ ...
  par.kappa_M).*par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.* ...
  par.eta_M.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1).* ...
  par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^((( ...
  -1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).* ...
  inp.K_B_R.^(((-1)+par.kappa_R).*par.kappa_R.^(-1)).*((-1)+ ...
  par.eta_R)+inp.L_R.*par.eta_R.*par.scale_L.*(inp.L_R.* ...
  par.scale_L).^((-1).*par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^( ...
  -1).*par.kappa_R)).^(((-1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu) ...
  .^(-1)).*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((-1) ...
  .*par.sigma.^(-1)).*((-1).*((-1)+par.mu).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^(((-1)+ ...
  par.rho).^(-1).*par.rho.*((-1)+par.sigma).*par.sigma.^(-1))+ ...
  par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^(((-1)+par.kappa_M).* ...
  par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.*par.eta_M.* ...
  par.scale_L.*(inp.L_M.*par.scale_L).^((-1).*par.kappa_M.^(-1))) ...
  .^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^(((-1)+par.nu).* ...
  par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+ ...
  par.kappa_R).*par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.* ...
  par.eta_R.*par.scale_L.*(inp.L_R.*par.scale_L).^((-1).* ...
  par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^((( ...
  -1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^((( ...
  -1)+par.sigma).*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+( ...
  -1).*par.alpha.*par.sigma));
        end        
        
        function out = Y_B_C(self,inp,varargin)

            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_B_C.^((-1).*par.kappa_C.^(-1)).*inp.K_S.^par.alpha.*par.A.*( ...
  (-1)+par.alpha).*((-1)+par.eta_C).*((-1)+par.lambda).*((-1)+ ...
  par.mu).*((-1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)) ...
  .*((-1)+par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1)).*(((-1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).* ...
  ((-1)+par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^((-1).*par.rho.^(-1)).*(inp.K_E.^(((-1)+ ...
  par.rho).*par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*((( ...
  -1).*inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^((par.rho+( ...
  -1).*par.sigma).*(par.sigma+(-1).*par.rho.*par.sigma).^(-1)).*(( ...
  -1).*((-1)+par.mu).*(inp.K_E.^(((-1)+par.rho).*par.rho.^(-1)).* ...
  par.lambda+(-1).*((-1)+par.lambda).*(((-1).*inp.K_B_C.^(((-1)+ ...
  par.kappa_C).*par.kappa_C.^(-1)).*((-1)+par.eta_C)+inp.L_C.* ...
  par.eta_C.*par.scale_L.*(inp.L_C.*par.scale_L).^((-1).* ...
  par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^(-1).*par.kappa_C)).^((( ...
  -1)+par.rho).*par.rho.^(-1))).^(((-1)+par.rho).^(-1).*par.rho.*(( ...
  -1)+par.sigma).*par.sigma.^(-1))+par.mu.*((par.upsilon.*(((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^(((-1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+ ...
  par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(((-1)+par.nu).* ...
  par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^(((-1)+par.sigma) ...
  .*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+(-1).*par.alpha.* ...
  par.sigma));
        end        
        
        function out = Y_E(self,inp,varargin)

            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_E.^((-1).*par.rho.^(-1)).*inp.K_S.^par.alpha.*par.A.*((-1)+ ...
  par.alpha).*par.lambda.*((-1)+par.mu).*(inp.K_E.^(((-1)+par.rho).* ...
  par.rho.^(-1)).*par.lambda+(-1).*((-1)+par.lambda).*(((-1).* ...
  inp.K_B_C.^(((-1)+par.kappa_C).*par.kappa_C.^(-1)).*((-1)+ ...
  par.eta_C)+inp.L_C.*par.eta_C.*par.scale_L.*(inp.L_C.* ...
  par.scale_L).^((-1).*par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^( ...
  -1).*par.kappa_C)).^(((-1)+par.rho).*par.rho.^(-1))).^((par.rho+( ...
  -1).*par.sigma).*(par.sigma+(-1).*par.rho.*par.sigma).^(-1)).*(( ...
  -1).*((-1)+par.mu).*(inp.K_E.^(((-1)+par.rho).*par.rho.^(-1)).* ...
  par.lambda+(-1).*((-1)+par.lambda).*(((-1).*inp.K_B_C.^(((-1)+ ...
  par.kappa_C).*par.kappa_C.^(-1)).*((-1)+par.eta_C)+inp.L_C.* ...
  par.eta_C.*par.scale_L.*(inp.L_C.*par.scale_L).^((-1).* ...
  par.kappa_C.^(-1))).^(((-1)+par.kappa_C).^(-1).*par.kappa_C)).^((( ...
  -1)+par.rho).*par.rho.^(-1))).^(((-1)+par.rho).^(-1).*par.rho.*(( ...
  -1)+par.sigma).*par.sigma.^(-1))+par.mu.*((par.upsilon.*(((-1).* ...
  inp.K_B_M.^(((-1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+ ...
  par.eta_M)+inp.L_M.*par.eta_M.*par.scale_L.*(inp.L_M.* ...
  par.scale_L).^((-1).*par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^( ...
  -1).*par.kappa_M)).^(((-1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+ ...
  par.upsilon).*(((-1).*inp.K_B_R.^(((-1)+par.kappa_R).* ...
  par.kappa_R.^(-1)).*((-1)+par.eta_R)+inp.L_R.*par.eta_R.* ...
  par.scale_L.*(inp.L_R.*par.scale_L).^((-1).*par.kappa_R.^(-1))) ...
  .^(((-1)+par.kappa_R).^(-1).*par.kappa_R)).^(((-1)+par.nu).* ...
  par.nu.^(-1))).^(((-1)+par.nu).^(-1).*par.nu)).^(((-1)+par.sigma) ...
  .*par.sigma.^(-1))).^(((-1)+par.sigma).^(-1).*(1+(-1).*par.alpha.* ...
  par.sigma));
        end

        function out = Y_S(self,inp,varargin)
            
            if self.parse_inputs
                par = self.parse_inputs_fun(varargin{:});
            else
                par = varargin{:};
            end
            
            out = inp.K_S.^((-1)+par.alpha).*par.A.*par.alpha.*((-1).*((-1)+par.mu) ...
  .*(inp.K_E.^(((-1)+par.rho).*par.rho.^(-1)).*par.lambda+(-1).*(( ...
  -1)+par.lambda).*(((-1).*inp.K_B_C.^(((-1)+par.kappa_C).* ...
  par.kappa_C.^(-1)).*((-1)+par.eta_C)+inp.L_C.*par.eta_C.* ...
  par.scale_L.*(inp.L_C.*par.scale_L).^((-1).*par.kappa_C.^(-1))) ...
  .^(((-1)+par.kappa_C).^(-1).*par.kappa_C)).^(((-1)+par.rho).* ...
  par.rho.^(-1))).^(((-1)+par.rho).^(-1).*par.rho.*((-1)+par.sigma) ...
  .*par.sigma.^(-1))+par.mu.*((par.upsilon.*(((-1).*inp.K_B_M.^((( ...
  -1)+par.kappa_M).*par.kappa_M.^(-1)).*((-1)+par.eta_M)+inp.L_M.* ...
  par.eta_M.*par.scale_L.*(inp.L_M.*par.scale_L).^((-1).* ...
  par.kappa_M.^(-1))).^(((-1)+par.kappa_M).^(-1).*par.kappa_M)).^((( ...
  -1)+par.nu).*par.nu.^(-1))+(-1).*((-1)+par.upsilon).*(((-1).* ...
  inp.K_B_R.^(((-1)+par.kappa_R).*par.kappa_R.^(-1)).*((-1)+ ...
  par.eta_R)+inp.L_R.*par.eta_R.*par.scale_L.*(inp.L_R.* ...
  par.scale_L).^((-1).*par.kappa_R.^(-1))).^(((-1)+par.kappa_R).^( ...
  -1).*par.kappa_R)).^(((-1)+par.nu).*par.nu.^(-1))).^(((-1)+par.nu) ...
  .^(-1).*par.nu)).^(((-1)+par.sigma).*par.sigma.^(-1))).^((-1).*(( ...
  -1)+par.alpha).*((-1)+par.sigma).^(-1).*par.sigma);
        end        
        
        function par = parse_inputs_fun(self, varargin)

            if self.KrusellRobotsEach_use_properties
                
                par.alpha = self.alpha;
                par.rho = self.rho;
                par.sigma = self.sigma;
                par.kappa_M = self.kappa_M;
                par.kappa_R = self.kappa_R;
                par.kappa_C = self.kappa_C;
                par.nu = self.nu;
                par.lambda = self.lambda;
                par.mu = self.mu;
                par.eta_M = self.eta_M;
                par.eta_R = self.eta_R;
                par.eta_C = self.eta_C;
                par.upsilon = self.upsilon;
                par.A = self.A;
                par.gamma = self.gamma;
                par.scale_L = self.scale_L;
                
            else
                
                p = inputParser;
                p.KeepUnmatched = true;
                
                p.addParameter('alpha', self.alpha, @(x)  x>0)
                p.addParameter('rho', self.rho, @(x) x>0)
                p.addParameter('sigma', self.sigma, @(x) x>0)
                p.addParameter('kappa_M', self.kappa_M, @(x) x>0)
                p.addParameter('kappa_R', self.kappa_R, @(x) x>0)
                p.addParameter('kappa_C', self.kappa_C, @(x) x>0)
                p.addParameter('nu', self.nu, @(x) x>0)
                p.addParameter('lambda', self.lambda, @(x) ((x>=0) && (x<=1)))
                p.addParameter('mu', self.mu, @(x) ((x>=0) && (x<=1)))
                p.addParameter('eta_M', self.eta_M, @(x) ((x>=0) && (x<=1)))
                p.addParameter('eta_R', self.eta_R, @(x) ((x>=0) && (x<=1)))
                p.addParameter('eta_C', self.eta_C, @(x) ((x>=0) && (x<=1)))
                p.addParameter('upsilon', self.upsilon, @(x) ((x>=0) && (x<=1)))
                p.addParameter('A', self.A, @(x) x>0)
                p.addParameter('gamma', self.gamma, @(x) x>0)
                p.addParameter('scale_L', self.scale_L, @(x) x>0)
                
                p.parse(varargin{:});

                tools.no_empty_fields(p.Results);
                par = p.Results;
            end
        end
        
    end    
    
end

